#pragma once 

#include "table.h"

#include <numbers>
#include <cmath>

namespace PS
{

template<int Dim>
class Point : public Table<1, scalar>
{
public:
    using value_type = scalar;
    
public:
    Point();

    virtual ~Point() = default;

    template<typename ... T>
    Point(T ... values);

    Point(const Point& other);

    Point(Point&& other) noexcept;

    Point& operator=(const Point& other);

    Point& operator=(Point&& other) noexcept;

    void set(scalar v1, scalar v2);

    void set(scalar v1, scalar v2, scalar v3);

    using TableBase<1, scalar>::reinit;

    scalar distance(const Point<Dim>& other) const;

    Point<Dim> rotate(scalar angle_in_angular) const;

    scalar length() const;

    Point<Dim> normalize() const;

    Point<Dim> projection(const Point<Dim>& other) const;

    scalar inner_product(const Point<Dim>& other) const;

public:
    Point<Dim> & operator+=(const Point<Dim> & other);

    Point<Dim> & operator-=(const Point<Dim> & other);

    Point<Dim> & operator*=(scalar value);

    Point<Dim> & operator/=(scalar value);

};

template<int Dim>
using Vector = Point<Dim>;

} // namespace PS


namespace PS
{

template<int Dim>
Point<Dim>::Point() : Table<1, scalar>(Dim)
{
    fill(scalar{0});
}

template<int Dim>
template<typename ... T>
Point<Dim>::Point(T ... values) : Table<1, scalar>(Dim)
{
    static_assert(Dim == sizeof...(T), "Wrong number of constructor arguments");
    set(values ...);
}


template<int Dim>
Point<Dim>::Point(const Point<Dim>& other)
{
    this->operator=(other);
}

template<int Dim>
Point<Dim>::Point(Point<Dim>&& other) noexcept
{
    this->operator=(std::move(other));
}

template<int Dim>
Point<Dim>& Point<Dim>::operator=(const Point<Dim>& other)
{
    if(&other == this) return *this;

    Table<1, scalar>::operator=(other);

    return *this;
}

template<int Dim>
Point<Dim>& Point<Dim>::operator=(Point<Dim>&& other) noexcept
{
    Table<1, scalar>::operator=(std::move(other));

    return *this;
}

template<>
inline void
Point<2>::set(scalar v1, scalar v2)
{
    AssertDimension(n_elements(), 2);

    (*this)[0] = v1;
    (*this)[1] = v2;
}

template<>
inline void
Point<3>::set(scalar v1, scalar v2, scalar v3)
{
    AssertDimension(n_elements(), 3);

    (*this)[0] = v1;
    (*this)[1] = v2;
    (*this)[2] = v3;
}

template <>
inline Point<2>
Point<2>::rotate(scalar angle_in_angular) const
{
    const auto angel_in_radian = std::numbers::pi * angle_in_angular / 180.0;

    return {((*this)[0] * std::cos(angel_in_radian) - (*this)[1] * std::sin(angel_in_radian)),
            ((*this)[0] * std::sin(angel_in_radian) + (*this)[1] * std::cos(angel_in_radian))};
}

template<int Dim>
scalar 
Point<Dim>::length() const
{
    scalar tmp = 0.0;
    for(int i = 0; i < Dim; ++i)
    {
        tmp += std::pow((*this)[i], 2);
    }
    
    return std::sqrt(tmp);
}


template<int Dim>
scalar 
Point<Dim>::distance(const Point<Dim>& other) const
{
    scalar tmp = 0.0;
    for(int i = 0; i < Dim; ++i)
    {
        tmp += std::pow((*this)[i] - other[i], 2);
    }
    
    return std::sqrt(tmp);
}

template<int Dim>
Point<Dim> 
Point<Dim>::normalize() const
{
    Point<Dim> rst = (*this);

    const auto len = length();

    for(int i = 0; i < Dim; ++i)
    {
        rst[i] /= len;
    }

    return rst;
}

template<int Dim>
Point<Dim> 
Point<Dim>::projection(const Point<Dim>& other) const
{
    const auto tmp = this->inner_product(other) / std::pow(other.length(), 2.0);

    return tmp * other;
}

template<int Dim>
scalar 
Point<Dim>::inner_product(const Point<Dim>& other) const
{
    scalar rst = 0.0;

    for(int i = 0; i < Dim; ++i)
    {
        rst += (*this)[i] * other[i];
    }

    return rst;
}


template<int Dim>
Point<Dim> & 
Point<Dim>::operator+=(const Point<Dim> & other)
{
    for(auto i = 0; i < this->n_elements(); ++i)
    {
        (*this)[i] += other[i];
    }

    return *this;
}

template<int Dim>
Point<Dim> & 
Point<Dim>::operator-=(const Point<Dim> & other)
{
    for(auto i = 0; i < this->n_elements(); ++i)
    {
        (*this)[i] -= other[i];
    }

    return *this;
}

template<int Dim>
Point<Dim> & 
Point<Dim>::operator*=(scalar value)
{
    for(auto i = 0; i < this->n_elements(); ++i)
    {
        (*this)[i] *= value;
    }

    return *this;
}    


template<int Dim>
Point<Dim> & 
Point<Dim>::operator/=(scalar value)
{

    Assert(std::isnormal(value), ExcDivideByZero());

    for(auto i = 0; i < this->n_elements(); ++i)
    {
        (*this)[i] /= value;
    }    

    return *this;
}

} // namespace PS


namespace PS
{

template<int Dim>
inline Point<Dim>
operator + (const Point<Dim>& left, const Point<Dim>& right)
{
    Point<Dim> rst(left);

    rst += right;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator + (const Point<Dim>& left, scalar right)
{
    Point<Dim> rst(left);

    rst += right;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator + (scalar left, const Point<Dim>& right)
{
    Point<Dim> rst(right);

    rst += left;

    return rst;
}


template<int Dim>
inline Point<Dim>
operator - (const Point<Dim>& left, const Point<Dim>& right)
{
    Point<Dim> rst(left);

    rst -= right;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator - (const Point<Dim>& left, scalar right)
{
    Point<Dim> rst(left);

    rst -= right;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator * (const Point<Dim>& left, scalar right)
{
    Point<Dim> rst(left);

    rst *= right;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator * (scalar left, const Point<Dim>& right)
{
    Point<Dim> rst(right);

    rst *= left;

    return rst;
}

template<int Dim>
inline Point<Dim>
operator / (const Point<Dim>& left, scalar right)
{
    Point<Dim> rst(left);

    rst /= right;

    return rst;
}

template <int Dim>
std::ostream &
operator<<(std::ostream &out, const Point<Dim> &point)
{
  out << '[';
  for (unsigned int i = 0; i < Dim; ++i)
    {
      out << point[i];
      if (i + 1 != Dim)
        out << ',';
    }
  out << ']';

  return out;
}

} // namespace PS
